Quasiequilibrium supersolid phase of a two-dimensional dipolar crystal 



O 

(N 

Ph. 
< 



I 

O 

o 



> 
o 

00 

o 

^' 

o 
o 



X 



I.L. Kurbakov/ Yu.E. Lozovik/ G.E. Astrakharchik,^ and J. Boronat^ 

^Institute of Spectroscopy, 142190 Troitsk, Moscow region, Russia 
^ Departament de Fisica i Enginyeria Nuclear, Campus Nord B4--B5, 
Universitat Politecnica de Catalunya, E-08034 Barcelona, Spam 
(Dated: April 7, 2010) 

We have studied the possible existence of a supersolid phase of a two-dimensional dipolar crystal 
using quantum Monte Carlo methods at zero temperature. Our results show that the commensurate 
solid is not a supersolid in the thermodynamic limit. The presence of vacancies or interstitials turns 
the solid into a supersolid phase even when a tiny fraction of them are present in a macroscopic 
system. The effective interaction between vacancies is repulsive making a quasiequilibrium dipolar 
supersolid possible. 

PACS numbers: 



Quantum systems with a dominant dipolar interaction 
have received permanent interest from the recent achieve- 
ment of a Bose-Einstein condensed state of chromium 
atoms with a large dipole moment [1]. The anisotropy of 
the dipole-dipole interaction leads to new exciting quan- 
tum phases that have been observed for the first time in 
fully quantum systems The experimental confirma- 
tion of these predicted phases and the collapses induced 
by the attractive part of the interaction can be even bet- 
ter realized if the permanent dipole moment of the par- 
ticles becomes larger. A promising system for this £oa\ 
is a stable gas of ultracold heteronuclear molecules [sf. 

If the quantum dipoles are confined in a two- 
dimensional (2D) plane, and all the dipole moments are 
perpendicular to the plane, the interaction is always re- 
pulsive and therefore the system is stable at any den- 
sity. Under such spatial and orientational restrictions 
one looses relevant features that can emerge when the at- 
tractive collapse is approached but, on the other side, the 
stability of the 2D geometry allows for the possible ob- 
servation of agas-solid quantum phase transition at high 
densities Q . A 2D environment is currently devised 
in the field of cold quantum gases by very anisotropic 
traps where the confinement in one direction is so tight 
that the transverse motion is frozen to zero-point oscil- 
lations [6] . Another physical system where this 2D setup 
is relevant is the one of indirect excitons composed by 
electrons and holes physically separated using two cou- 
pled quantum wells 0, |1] . If the distance between the 
electron and hole layers is significantly smaller than the 
electron-electron and hole-hole distances the resulting ex- 
citons can be modeled as composite bosons with a dipole- 
dipole interaction [13] ■ 

Recent quantum Monte Carlo calculations at zero Q 
and finite temperature [5] have shown that a 2D homo- 
geneous phase of dipoles experiments a gas-solid phase 
transition when the density increases. The equation of 
state of this solid phase, which forms a triangular lattice, 
as well as its main structure properties are already re- 
ported in these previous works. However, relevant ques- 



tions such as the possible superfluid signal and/or con- 
densate fraction of the zero-temperature dipolar crystal 
were not addressed so far. In fact, there is at present a 
renewed interest in the search of supersolid phases where 
ofF-diagonal long-range order and spatial solid order are 
simultaneously present flU . Supersolids appear as a new 
intriguing state of matter that was predicted long-time 
ago and only recently observed in torsional oscillator ex- 
periments with solid ^He [l2j . 

Within the framework of Bose-Hubbard Hamiltonians 
supersolid phases of dipolar lattice bosons have already 
been predicted. Danshita and Sa de Melo iden- 
tify exotic phases as checkerboard and striped supersolid 
phases by including in the model Hamiltonian the attrac- 
tive part of the dipole-dipole interaction and Trefzger et 
al. [14] find a pair-supersolid phase in a bilayer configu- 
ration. The emergence of supersolid states when dipolar 
bosons are confined in two-dimensional optical lattices is 
probably favored by the free tuning of the localization 
strength of the external lattice potential included in the 
model Hamiltonian. A different concern is the possible 
formation of a supersolid phase in a continuum system 
where a solid is formed at high density without the pres- 
ence of any external localization potential. In this work, 
we present the first study of supersolidity in 2D dipolar 
bosons at zero temperature using quantum Monte Carlo 
methods that rely merely on the microscopic Hamilto- 



The triangular crystal phase of dipolar bosons is stud- 
ied by means of the diffusion Monte Carlo (DMC) 
method that it is nowadays a standard tool for achieving 
exact solutions of many-boson systems at zero tempera- 
ture relying on a stochastic approach [3l • The Hamilto- 
nian describing the 2D system of N dipoles is 
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m being the mass and V{r) the dipole-dipole potential 
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The constant Cdd depends on the nature of the dipole- 
dipole interaction and increases proportionally to the 
square of the individual dipole moment. As in previ- 
ous work, we define characteristic units of length tq = 
™C'dd/(47r?i^) and energy £q = /{mr'^) in such a way 
that the properties of the system are governed by a di- 
mensionless density nrg, with n the particle density. In 
order to fix the symmetry of the system and reduce the 
statistical variance DMC introduces a trial wave function 
^ that is used for importance sampling. In a previous 
work [ij, the crystal phase of dipoles was studied using 
for ^ a non symmetric model (Nosanow-Jastrow (NJ)) 
since the focus was the determination of the equation of 
state and the phase transition point, issues in which im- 
plicit symmetrization is much less relevant. Obviously, 
the NJ trial wave function can not be used in the present 
study since our goal is the determination of superfluid 
signals in the solid and that is only possible assuming 
particle indistinguishability. To this end, in the present 
work we use a symmetric model 
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that was first introduced in the study of solid ^He at zero 
temperature pJi]. In Eq. ([3|), R ~ {ri, . . . , rjv}, /(r) is a 
two-body Jastrow correlation factor chosen as in Ref. [3] , 
gifii) = exp[— Q!(ri — r/)^], and A'cr is the number of lat- 
tice sites of the triangular crystal structure. This model 
wave function ([3]) makes compatible the spatial solid or- 
der and the symmetry under the interchange of particles 
avoiding the numerically unworkable use of permanents 
on top of the NJ wave function. 

Coherence phenomena in the dipolar solid have been 
studied by calculating the one-body density matrix pi (r) 
and the superfluid fraction risjn. The function p\(r) ap- 
proaches a constant at long distances, which is the con- 
densate fraction Nq/N — lim^^oo Pi{T)/n, if off-diagonal 
long range order exists in the system. In DMC, the func- 
tion pi{r) can not be calculated using a pure estimator 
and therefore some bias induced by the trial wave func- 
tion 4' remains. To reduce this bias as far as possible the 
variational parameters entering in ^' have been optimized 
in such a way that the variational and DMC (mixed) es- 
timations of Nq/N are coincident within their statistical 
errors. On the other hand, the superfluid density is com- 
puted by extending the winding-number technique, used 
in path integral Monte Carlo (PIMC) simulations at fi- 
nite temperature, to zero temperature [13] • Explicitly, 
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FIG. 1: (Color online) (a) One-body density matrix of the 
perfect 2D crystal at a density nrg = 290 and as a function 
of the number of particles used in the simulation; the solid 
line corresponds to a nonsymmetric trial wave function, (b) 
The function aDs{T), Eq. 21 at the same density and as a 
function of A'^; the thin straight line corresponds to Us/n — 1 
and the solid line to the nonsymmetric case, Ua/n — 0. 



where a = N/{ADq) with Dq = h'^/{2m), and Dsir) = 
((Rcm(''') — R.cm(0))^), with RcM the center of mass of 
the particles and r the imaginary time. Differently from 
the estimation of pi{r), the measure of the superfluid 
density (HI) is unbiased (pure estimator). 

DMC results for the perfect 2D triangular solid are re- 
ported in Fig. [T] In all the simulations, carried out with 
different number of particles N, the one-body density 
matrix shows a plateau at long distances and therefore 
a finite condensate fraction. However, Ng/N decreases 
significantly with N making the condensate fraction van- 
ishingly small in the thermodynamic limit N ^ oo. If 
the calculation is carried out with a nonsymmetric trial 
wave function (Nosanow-Jastrow model), pi{r) does not 
show off-diagonal long range order for any value of TV (see 
Fig. [T|) . We show in the same figure results for the super- 
fiuid density, plotting the function aDgir) (Eq. [3]) as a 
function of the imaginary time r; the slope of this func- 
tion is directly Us/n. As one can see, the slope becomes 
zero within our numerical resolution for values iV > 30 
pointing to the absence of supersolidity in the perfect 
crystal in the thermodynamic limit. The lack of super- 
solid signatures in the commensurate solid is observed at 
any density, starting on the melting one T^m^'o — 290(30) 
shown in Fig. [T] 

The presence of defects or imperfections in a crystal 
has been suggested as a plausible explanation of the su- 
perolid signals observed experimentally in torsional oscil- 
lator measurements of solid ^He. Whereas there are still 
open discussions about the existence or not of vacancies 
in the ground-state of bulk solid ^He [ll| , it seems very 
reasonable to think on the presence of vacancies or in- 
terstitials in a 2D crystal of dipoles. Indeed, indirect 
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FIG. 2: (Color online) Condensate fraction (a), superfluid 
fraction (c), and height of the main divergent peak of 
AS(fciatt) with respect to the background at the reciprocal 
lattice vector fciatt (e), for a crystal with vacancies or inter- 
stitials (A^ ^ Ncr)- Density dependence of the condensate 
fraction (b), superfluid fraction (d), and S{k) peak (f) for a 
solid with one vacancy. Solid and open points stand for the 
solid and gas phases, respectively. 



excitons or trapped atoms or niolecules with high dipole 
moments can be easily produced with a fraction of de- 
fects. We have collected in Fig.[2]DMC results for a dipo- 
lar solid with a finite fraction of vacancies or interstitials. 
Both the condensate fraction and superfluid density are 
very sensitive to the presence of defects: for any concen- 
tration of vacancies or interstitials a finite value for Us/n 
and Nq/N is observed. When the fraction of vacancies 
is ^ 6% and nr^ = 290, the supersolid melts: Nq/N 
equals its value in the gas phase, ng/n = 1, and the peak 
of S{k) in the reciprocal lattice vector disappears. The 
crystal also melts due to interstitials at a slightly higher 
concentration, ~ 10%. 

In Fig. [21 we show the density dependence of Nq/N, 
Hs/n, and height of the main peak of S{k) for the par- 
ticular case of one vacancy in a solid with iVcr = 30. The 
condensate fraction becomes vanishingly small at high 
densities and remains always a factor of 3-5 smaller than 
its value in the gas phase, as also shown in the figure for 
comparison. The superfiuid density fraction is ~ 50% at 
melting of the commensurate crystal and decreases with 
n but much more slowly than Nq/N. On the other hand, 
the height of the S{k) peak increases with density as ex- 
pected, and it increases with N for a fixed n as it must 
be in a solid structure (not shown in the figure). At den- 
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FIG. 3: (Color online) (a) Snapshot of a typical spatial con- 
figuration when vacancies are present. This case corresponds 
to A'cr = 90 and A'' — 86; crosses are the lattice points, solid 
circles are particles, and open circles, vacancies, (b) Vacancy- 
vacancy pair distribution function gvv{r) for the same A'cr and 
A' values as in (a); L = \jN/n is the size of the simulation 
box. 



sity nr\ — 230 the supersolid completely melts: Nq/N 
becomes equal to its value in the gas, Us/n = 1, and the 
divergent peaks in Sik) disappear. 

As we commented before, the condensate fraction 
shows a significant dependence with N and therefore the 
estimation of the thermodynamic limit when vacancies 
are present is fundamental. For this purpose, we per- 
formed a study of the iV-dependence of Nq/N for va- 
cancy fractions 0.018 < iVvac/iVcr < 0.042. The DMC 
results obtained show a \/N decrease with the number 
of Bose-condensed particles per vacancy, Nq/N^hc but 
with a finite value in the thermodynamic limit iV — )■ oo 
of A^o/iVvac = 0.050(8). The number of superfluid parti- 
cles per vacancy, which is weakly dependent on TV, also 
remains finite in this limit. Also, the height of the nar- 
row peak in S{k) remains finite and proportional to N. 
Therefore, vacancy-induced superfluidity coexists with 
spatial solid order, i.e., a supersolid phase can exist. 

A relevant concern about the stability of a small frac- 
tion of vacancies in the solid is the nature of their mu- 
tual interaction. Several microscopic estimations in solid 
■*He show that two vacancies tend to form a weakly- 



bound state because their interaction is attractive [18- 
[ioj . Therefore, it has been argued that vacancies would 
form a cluster inside the crystal that eventually can evap- 
orate producing a collapse of the crystal. In order to 
characterize the local structure of vacancies in a 2D dipo- 
lar crystal we have sampled the vacancy-vacancy two- 
body distribution function gwir). As vacancies are not 
real particles and our simulation works in a configura- 
tion space of particle coordinates one has to define what 
a vacancy position is for a given snapshot of the system. 
In our procedure, we have always identified a vacancy 
with one of the sites of the perfect triangular lattice in 
which unambiguously none of the particles is around it 
within a cutoff radius that is close to the value of the 
lattice constant. Along the evolution in imaginary time, 
there are configurations in which we can not identify the 
vacancy sites due to intrinsic fluctuations; in these cases 
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FIG. 4: (Color online) Energy per particle of the 2D dipolar 
solid as a function of the number of vacancies or interstitials, 
at density nrg — 290 and Na — 90, and normalized to the 
energy per particle of the commensurate crystal. 



we do not accumulate statistics for gvv{r). In Fig. [3l we 
show a snapshot of the system where vacancy sites are 
identified according to our definition. In the same fig- 
ure, the vacancy- vacancy correlation function gwir) is 
shown for a setup composed by TVcr — 90 and four va- 
cancies. The radial function gwif) is normalized at each 
distance dividing what accumulated in each bin by the 
corresponding output obtained in a merely random dis- 
tribution. As shown in Fig. |3l vacancies repel at short 
distances and this relevant feature is not only observed 
for this particular set of parameters but settled at other 
densities and vacancy fractions (below the threshold for 
melting). The repulsive interaction between vacancies in 
the 2D dipolar solid is probably due to the monotonously 
repulsive interaction between aligned dipoles, that makes 
configurations be more stable when vacancies spread in 
the system in order to effectively reduce the dipolar den- 
sity. This is contrary to the vacancy- vacancy attraction 
observed in solid *He simulations in which the van der 
Waals attraction at long distances can explain the differ- 
ence with the present results. 

The ground state of the dipolar solid is a commensu- 
rate phase, i.e., without vacancies and/or interstitials. In 
other words, an activation energy is needed to create a 
vacancy or interstitial. We have estimated the activation 
energy to create one or more vacancies or interstitials us- 
ing the DMC method. In Fig. lU we show the energies 
per particle of a solid with defects normalized to the en- 
ergy per particle of the commensurate solid and having 
rescaled the simulation box size to work at fixed density. 
Our results support the metastability of the solid with 
defects with respect to the commensurate solid. The ac- 
tivation energy for the creation of a vacancy is higher 
than the corresponding one for an interstitial and this 
is maintained when the number of defects increases: the 
slopes of the two cases are rather different (see Fig. |4]) . 
It is worth noticing that this behavior is opposite to the 
one observed in solid "^He where the activation energy for 
an interstitial is significantly larger than for a vacancy. 
Defining the activation energy in the standard way 21| . 
we get for one vacancy = 2150(150) and for one in- 



terstitial Ei = 990(50), both at fixed density nr^ = 290. 
These activation energies are significantly larger than the 
Berezinskii-Kosterlitz-Thouless temperature T = 380 [H] , 
making thermal activation of defects in a dipolar crystal 
not possible. 

Summarizing, we have studied the possible emergence 
of bosonic coherence phenomena in a two-dimensional 
crystal of dipoles by calculating the condensate fraction 
and superfluid density using accurate quantum Monte 
Carlo methods. To this end, we have used for the first 
time in this system a trial wave function for importance 
sampling with both boson symmetry and solid order. 
Our DMC results show that the commensurate solid is 
not a supersolid since both ris/n and Nq/N become zero 
in the thermodynamic limit within our numerical resolu- 
tion. The introduction of defects, in the form of vacancies 
or interstitials, produces a dramatic effect on both quan- 
tities, even with tiny concentrations. A quasiequilibrium 
solid with vacancies or interstitials is proven to be su- 
persolid within a predicted fraction of defects. If this 
percentage is further increased the supersolid melts. The 
effective vacancy- vacancy interaction is repulsive at short 
distances, a feature that is opposite to the one of solid 
^He and that can help to stabilize the dipolar supersolid 
phase. The recrystallization to the ground state is expo- 
nentially suppressed by the tunelling barrier which stabi- 
lizes a dipolar supersolid with a small fraction of defects. 
Possible experimental realizations of a quasiequilibrium 
dipolar supersolid with defects include i) harmonically 
trapped dipolar molecules [131 or atoms ^23] spatially 
localized with an optical lattice and ii) Wannier-Mott 
2D dipolar excitons in single or coupled semiconductor 
quantum wells in electric and magnetic fields which are 
perpendicular to the quantum wells plane [j^ . In the 
latter case, the finite excitation lifetime caused by their 
optical recombination gives rise to a continuous addition 
of vacancies into the system, resulting in a macroscopic 
supersolid at low temperatures. 
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